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Abstract 


Managed aquifer recharge (MAR) is a water management strategy that 
uses aquifers for the seasonal or inter-annual storage of exceeding 
surface water. Systems that implement artificial recharge use infiltration 
ponds or shallow or deep infiltration wells. While articifial recharge 
through infiltration ponds requires large areas and infiltration through 


deep wells demands high investment, infiltration through drywells, which 
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inject water to the vadose zone, is often a preferred option because of its 
low cost and relatively easy implementation. We present results of 
detailed numerical simulations to assess the operation of a single drywell 
to implement an artificial recharge system. We demonstrate that 
depending upon the distribution and property of the geological materials 
that compose the subsurface, the operation of the well can produce quite 
different results in terms of effective recharge. We expect that these 
findings will help design real MAR systems and provide supporting 
information for decision-makers in charge of approving and financing such 


water management systems. 


Keywords: Unsaturated flow, numerical simulations, artificial recharge, 
aquifer recharge, water management, aquifer storage and recovery 
(ASR). 


Resumen 


La recarga artificial (managed aquifer recharge, MAR) es una estrategia 
de gestión de los recursos hídricos que utiliza a los acuíferos para el 
almacenamiento temporal de excedentes de agua superficial. Los 
sistemas que implementan la recarga artificial de acuíferos usan piscinas, 
o pozos profundos o someros para infiltrar agua. Mientras la recarga 
artificial a través de piscinas requiere grandes extensiones de terreno y 
la recarga mediante pozos profundos demanda altos costos de inversión, 
la infiltración mediante pozos secos —que inyectan agua a la zona no 
saturada— es usualmente preferida como alternativa de recarga debido a 


su bajo costo y simple implementación. En este artículo presentamos 
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resultados de simulaciones numéricas detalladas para evaluar la 
operación de un sistema de recarga artificial compuesto de un solo pozo 
seco. Demostramos que dependiendo de la distribución y de las 
propiedades de los sedimentos que componen el subsuelo, la operación 
del pozo puede producir diferentes resultados en términos de recarga 
efectiva. Esperamos que estos resultados sean útiles para el diseño de 
sistemas de recarga artificial y sirvan de soporte para tomadores de 
decisión a cargo de aprobar y financiar sistemas similares de gestión de 


agua. 


Palabras clave: flujo no saturado, simulaciones numéricas, recarga 
artificial, recarga de acuíferos, gestión del agua, almacenamiento de agua 


en acuíferos y recuperación. 
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Introduction 


Economic development and population growth combined with the effects 
of climate change have resulted in increasing stress on groundwater 
systems, particularly in arid and semi-arid regions of the world (Lall, 
Josset, € Russo, 2020). The implementation of artificial or managed 
aquifer recharge (MAR) systems has gained interest in recent decades as 


a water management strategy (Kimrey, 1989; Sharma, Sharma, Mehta, 
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8 Marwaha, 2000; Gale, Neumanmn, Calow, € Moench, 2002; DGA, 2013; 
Dillon et a/., 2019; CNR, 2020; Ebrahim, Lautze, € Villholth, 2020; 
Algahtani, Sale, Ronayne, 8 Hemenway, 2021). The basic principle of 
such systems is using the storage space available in the subsurface to 
store water between periods of exceedance and shortage (Bouwer, 1996; 
Bouwer, 2002). The successful implementation of artificial recharge 
projects requires (Bouwer, 1996; Bouwer, 2002): water that can be 
stored during at least a specific period, storage volume in the subsurface, 
and infiltration capacity to permit injecting a significant amount of water 
with small effort. Often, the implementation of MAR systems aim at using 
the recharged water for an environmental purpose, e.g. the recovery of 
depleted aquifers or springs. Alternatively, recharged water is recovered 
for a different use, in which case they are known as aquifer storage and 
recovery (ASR) systems (Pyne, 2017; Algahtani et al., 2021). 
Additionally, artificial recharge serves to recover overexploited aquifers, 
prevent subsidence underneath urban areas, prevent or control saline 
intrusion in coastal aquifers, dispose of excess water due to dewatering 
during construction projects and store treated wastewater (NRC, 1994; 
Wilson et al., 1995; Dillon et a/., 2019; van Lopik, Hartog, € Schotting, 
2020). 


There are several technical options to enhance recharge that can be 
broadly classified into three categories: surface systems such as ponds or 
trenches, shallow wells that penetrate only into the vadose zone (also 
known as drywells), and deep wells that extend below the water table and 
reach the aquifer as shown in Figure 1 (Bouwer, 2002). Infiltration ponds 


or trenches have been implemented in many locations worldwide due to 
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their simplicity and low cost. However, these systems require extensive 
land area to attain significant infiltration rates, thus, they are difficult to 
implement in urban areas where land is relatively scarce and expensive. 
Moreover, the adequacy of infiltration ponds to produce effective recharge 
may be relatively low in locations where the water table is deep because 
of the occurrence of perched aquifers and the long time that may take 
water to reach the water table (Flint, Ellett, Christensen, € Martin, 2012). 
Drywell systems are also relatively inexpensive to build, since in most 
cases, they are dug by hand or using inexpensive drilling equipment or, 
are installed using direct push (Hándel, Liu, Dietrich, Liedl, € Butler, 
2014). Since the infiltration rate depends on the height of the water 
column inside the well, they can provide good infiltration capacity 
compared to shallow systems implemented with ponds (Bandeen, 1988; 
Edwards, Harter, Fogg, Washburn, € Hamad, 2016; Justino, Failache, 8 
Barbassa, 2021). However, the infiltration capacity of shallow wells may 
degrade due to the accumulation of fine sediments carried by the 
infiltrated water and/or the development of a biological layer around the 
well casing (bioclogging). This clogging creates an operational problem 
due to the impossibility of applying some re-development techniques such 
as periodic pumping (Bouwer, 2002; Glass, Simúnek, 8. Stefan, 2020). 
Finally, the implementation of deep well systems may require 
considerably more technical skill and the use of more expensive drilling 
equipment, which usually results in much higher implementation and 
operational costs (Bouwer 1996; Bouwer, 2002; Gale et al., 2002; Dillon 
et al., 2019). 
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Figure 1. Schematic of artificial recharge systems: trenches or ponds 
(A), infiltration wells that infiltrate into the vadose zone known as 
drywells (B), and deep wells that infiltrate directly into aquifers (C). 
Materials that composed the subsurface are classified into 
hydrostratigraphic units (HU) that exhibit certain uniformity. 


Infiltration ponds and shallow infiltration wells inject water into the 
unsaturated zone, hence they allow for natural treatment of the infiltrated 
water due to the filtration capacity of soils and the occurrence of chemical 
and biological processes, such as denitrification (Sasidharan, Bradford, 
Simúnek, 8: Kraemer, 2020; Gorski, Dailey, Fisher, Schrad, 8 Saltikov, 
2020; Hágg et al., 2021). Therefore, this type of system is preferred in 


some cases due to the potential enhancement of the recharged water 
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quality (Bouwer, 1996; Edwards et al., 2016; Dillon et a/., 2019). For 
example, they are used to implement soil aquifer treatment systems, SAT 
(e.g. Wilson et a/., 1995). This enhancement can be a major consideration 
for implementing MAR systems in countries or areas where regulation is 


strict concerning the recharged water quality. 


The infiltration rate into shallow unsaturated soils (vadose zone) 
depends on the saturated hydraulic conductivity of the sediments (K), the 
hydraulic head at the infiltration point (water depth and elevation), and 
the negative pressure head (h) at the tip of the saturation front, ¡.e. the 
volume of soil within which moisture content increases due to infiltration 
(Bouwer, 1996; Bouwer, 2002). For wells, it also depends on the effective 
radius; wells with a larger diameter can infiltrate more water. Water 
moves preferentially downward from the infiltration point through a 
homogeneous material until it reaches a material of lower permeability 
that may act as a flow barrier inducing a lateral or horizontal movement 
of water (Hándel et a/., 2014; Maples, Fogg, € Maxwell, 2019; Wu et al., 
2021). The net result is a much slower overall vertical water movement 
across the soil column and, hence, less effective recharge. If the source 
of the infiltration stops, water accumulated within the initially unsaturated 
sediments drains under the action of gravity. This drainage process may 
take a long time to complete after the infiltration source stops providing 
water (Bouwer, 1996; Algahtani et a/., 2021). The reduction in the vertical 
infiltration capacity produced by lower permeability units results in a 
leaking system composed of saturated materials overlying higher 
permeability units that are partially saturated (Flint et a/., 2012). Similar 


partially saturated leaking flow configurations occur in natural systems, 
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e.g. beneath river systems with low permeability sediments covering the 


streambed (Brunner, Cook, 8 Simmons, 2009). 


Infiltrated water may take a long time, of the order of years to 
decades, to reach the water table and contribute to recharge in presence 
of aquitards or hydrostratigraphic units with reduced K. This may render 
the recovery of the infiltrated water difficult (Bouwer, 1996; Flint et al., 
2012; Sasidharan, Bradford, Simúnek, 8: Kraemer, 2019; Sasidharan et 
al., 2020; Algahtani et a/., 2021). In addition, the recovery can be 
problematic due to difficulties in pumping groundwater from perched 
aquifers that might develop around the infiltration point (Bouwer, 2002). 
Such a long transit time can reduce the interest in this type of project by 
investors that seek returns in a limited time frame (Gale et a/., 2002). 
Therefore, selecting one or another system must weight several technical 
and economic factors (Table 1). When new systems are implemented 
using ponds or short infiltrations wells that inject into the vadose zone 
due to lower costs, lower technical barriers and potential advantages for 
improving the quality of the infiltrated water; the potential delay of the 
benefits due to the delayed recharge should also be accounted for (Gorski 
et al., 2020; Hágg et al., 2021). 
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Table 1. Factors considered for the design and operation of artificial 
recharge systems. Alternative methods correspond to infiltration pond 
(A), dry (B) and deep (C) infiltration wells. 


E Low Easy Easy Natural Water 
Description : . . ] 
cost [construction operation |filtration | recovery 
Infiltration pond X X X X 
Drywell X X X 
Deep well X 


Source: Bouwer (2002), Dillon et a/. (2019), Edwards et a/. (2016), Gale 
et al. (2002), Hágg et al. (2021), NRC (1994), Sharma et al. (2000), 
Wilson et al. (1995). 


Drywells have been analyzed in detail in several studies due to their 
advantages in terms of: construction cost, use of space in urban or well- 
developed agricultural areas, enhanced infiltration capacity in comparison 
to ponds or trenches, and their ability to infiltrate water into the vadose 
zone (Bandeen, 1988; Hándel et a/., 2014; Edwards et al. 2016; Liang, 
Zhan, € Zhang, 2018; Sasidharan et a/., 2019; Sasidharan et al., 2020; 
Fuentes, Chávez, Quevedo, Trejo-Alonso, 8 Fuentes, 2020; Glass et al., 
2020; van Lopik et a/., 2020; Justino et a/., 2021). Bandeen (1988) used 
numerical simulations to study the operation of a single drywell for the 
infiltration of excess rainfall considering a layered aquifer and a relatively 
shallow water table (30 meters deep). He simulated infiltration for a few 
hours after the beginning of a single rainfall event of 1 hr duration, 


considering only a small area around the well (< 30 meters). Hándel et 
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al. (2014) used numerical simulations to evaluate the potential 
performance of small diameter shallow infiltration wells installed with 
direct push. They demonstrated that drywells are cost-effective versus 
infiltration ponds because of the higher infiltration capacity that they may 
attain. They concluded that the presence of low permeability units 
significantly reduces the infiltration capacity. They considered only 
systems with a shallow water table so that infiltrated water always 
reached the aquifer. Edwards et al. (2016) review the literature related 
to drywells and comment on their advantages and limitations. They 
conclude that drywells are effective to enhance natural recharge to 
aquifers, but that their performance depends on well maintenance, 
quantity and quality of available water for recharge and local subsurface 
conditions. Liang et al. (2018) derive analytical solutions to model water 
infiltration produced by drywells and compare them with numerical 
simulations considering uniform materials. Sasidharan et al. (2020) 
performed numerical simulations of the operation of a drywell considering 
homogeneous and heterogeneous materials in a 2D axisymmetric domain 
delimited by a distance of 48 meters from the injection well. They 
evaluated ten realizations of heterogeneous materials generated with a 
geostatistical model that considered variance and spatial correlation in 
the horizontal and vertical directions and a water table located 20 meters 
below the injection point. They concluded that the arrival time and 
location of the recharge water can be substantially different depending on 
the heterogeneity of the materials. Glass et a/. (2020) proposes a model 
to incorporate variations in hydraulic conductivity due to clogging during 


the operation of MAR systems. They used the model to simulate two 
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laboratory experiments and evaluate several hypothetical recharge 
systems. They concluded that clogging can produce significant reductions 
in infiltration rates and large increases in pressure head that may create 


problems for the operation of real systems. 


Previous studies demonstrated that the performance of drywells in 
sites where the shallow water table can be controlled, at least locally 
around the injection well, by the distribution and characteristics of the 
sediments. We aim at extending the previous analyzes to a larger scale 
and to sites where the water table is deeper. We perform detailed 
numerical simulations to analyze the performance of a single infiltration 
well for different scenarios defined by the stratigraphy of the sediments 
located between the infiltration point and the water table. The level of 
detail of the simulations and the unsaturated flow modeling allow us to 
make a realistic evaluation of this kind of recharge system. Unlike 
previous studies, we simulate a large volume around the infiltration well 
that extends a few hundred meters from the well and considers a deep 
water table. We simulate the operation of a relatively deep drywell (> 50 
m) to represent real conditions of wells specifically constructed for 
infiltration in aquifers with a deep water table (> 50 m). These wells can 
take advantage of a longer water column to attain high well head pressure 
and, hence, high infiltration rates. We use the results of the numerical 
simulations to estimate infiltration rates and transit times to reach the 
water table. The main objective of this article is assessing potential 
advantages and limitations of recharge systems based on drywells 
considering conditions similar to a real-world setting. We expect that the 


findings can be useful for the design and implementation of similar 
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systems and, hence, will contribute to the further development of artificial 


recharge projects. 


Numerical simulations 


Setup 


The conceptual model for the simulation setup considered field conditions 
found in Santiago, Chile, which corresponds to a semi-arid climate 
(average rainfall < 300 mm/year). There, the water table is relatively 
deep and sediments are a mixture of alluvial sediments that originate 
from the Andes and fluvial deposits near rivers. Figure 2 shows the 
stratigraphic sequence observed in a more than 80 meters deep gravel 


pit, which is a good example of the sediments found in the area. 
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Figure 2. Stratigraphic sequence observed in a gravel pit near the 
Maipo River, Santiago, Chile, which served as basis for the 
conceptualization of the numerical simulations. The pit is approximately 


80-90 m meters deep. 
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We consider a site where the water table is more than 100 m deep 
and sediments are composed of a mixture of clays, sandy and gravely 
materials that form a relatively uniform stratigraphic column with not 
easily recognizable hydrogeological units (HU). For a similar geological 
setting, Sudicky (1986) reports changes in K measured through 
permeameter tests of 5 cm undisturbed soil samples of more than one 
order of magnitude within 2 meters of a shallow sandy aquifer located in 
Borden, Canada (Figure 5 in Sudicky, 1986). The sediments of the Borden 
aquifer were described as: “primarily horizontal, discontinuous lenses of 
medium-grained, fine-grained, and silty fine-grained sand”, based on a 
very detailed laboratory characterization of thousands of soil samples 
(Sudicky, 1986). This characterization of the sediments was performed as 
part of a large-scale field experiment that had as main purpose to 
understand the role of natural geological heterogeneity for solute 
transport in groundwater (Sudicky €: TIllman, 2011). The level of detail and 
sophistication used for the characterization of the Borden site sediments 
are well beyond the one routinely applied in hydrogeological projects. 
Identifying different hydrogeological units in similar geological settings is 
usually difficult to achieve based on borehole logging alone, which is the 


standard method used in most practical applications. 


We want to evaluate the potential impact that mild heterogeneity, 
such as the one reported for the Borden site, could produce for the 
practical implementation of an artificial recharge system based on shallow 
infiltration wells that do not reach the water table. Hence, we divide the 


model domain into three main horizontal hydrogeological units defined by 
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different hydraulic conductivity values (Table 2), which allows us to assess 
the potential impact of the presence of HUs with different hydraulic 
properties. This conceptualization is akin to the multilayered aquifer 
system used in many hydrogeological studies (e.g. Matheron € De 
Marsily, 1980; Fuentes et a/., 2020). The hydraulic conductivity values 
assigned to the different HUs of the model have relatively minor changes 
between different HUs, with a maximum ratio between the most and least 
permeable units of less than 100. For all HUs, we assume an anisotropy 
ratio K,/K, = 10, hence resistance to horizontal flow is less than to vertical 
flow. We evaluate six different scenarios based on the combination of K 
values for each one of the three HUs: One that considers the presence of 
only one homogenous material (S1) with a single value of K and five 
heterogeneous (S2 to S6) that consider different values of K for each HU 
(Table 2). For S1, we assume a value of K equal to the effective value for 
flow perpendicular (vertical) to the layered units of scenario S2. The 
unsaturated flow properties (see Appendix for details), ¡.e. water 
retention curve and relative permeability, were set equal to the 
corresponding to a loamy sand soil using the van Genuchten model (van 
Genuchten € Nielsen, 1985), and were considered equal for all three 
hydrogeological units. Hence, we only tested the influence of considering 
different values of saturated hydraulic conductivity. Considering different 
retention and relative permeability curves would result in relative minor 


differences in the results discussed below. 
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Table 2. Scenarios simulated assuming three hydrogeological units 
(HU): upper and lower sand-gravel units and middle sand with fine 
sediments unit, Scenarios S5 and S6 are similar to S4 except by the 
extension of the middle low permeability unit, which extends to only 300 
(S5) and 550 (S6) meters from the injection well. Porosity for all layers 
was set equal to 0.41, while initial moisture was set equal to the 


residual moisture content of a loamy sand, 0.057. 


a y Top Depth Bottom Kh Kv 
Scenario Unit 
(m) Depth (m) | (m/d) (m/d) 

si Single 0 230 43.30 4,33 
Upper 0 90 18.20 1.82 

S2 
Lower 90 230 60.70 6.07 
Upper 0 80 18.20 1.82 
S3 Fine 80 90 1.10 0.11 
Lower 90 230 60.70 6.07 
Upper 0 80 18.20 1.82 
S4 Fine 80 90 0.11 0.01 
Lower 90 230 60.70 6.07 


We simulate the operation of a single injection well using VS2DT 
(Lappala, Healy, 8 Weeks, 1987), which simulates unsaturated flow. 
VS2DT solves Richards” equation formulation based on pressure head 
using a finite difference approximation with a fully implicit backward time 
approximation. The non-linear partial differential equation is linearized 


using a Newton-Raphson scheme, and relative permeability and hydraulic 
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conductivity are evaluated using upstream weighting and harmonic 


average, respectively (Lappala et al., 1987). 


We used a two-dimensional grid assuming axial symmetry so that 
the coordinates of the grid correspond to the distance from the injection 
well, r, and the elevation, z. Figure 3 shows a schematic of the model 
setup, including boundary conditions. The numerical grid is equivalent to 
a cylindrical domain that extends 600 m from the injection well located at 
r = 0 m. To simplify the analysis, we consider that the origin of the vertical 
coordinate system (z = O m) coincides with the position of the ground 
surface and that it points downward. Vertically, the domain includes two 
zones: a fully saturated of height Hsat = 90 m and an unsaturated height 
Hunsat = 140 m above the water table. Including a saturated zone below 
the water table explicitly simulates a more realistic interaction between 
vertical infiltration and groundwater than assuming a free drainage 
boundary condition as used in other studies (e.g. Glass et al., 2020; 
Sasidharam, 2019, 2020). In particular, it allows evaluating changes in 
piezometric levels that could be monitored in real recharge systems to 


evaluate the magnitude and time evolution of the effective recharge. 
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Figure 3. Schematic (not drawn to scale) of simulated artificial recharge 
system. Injection well (left) simulated as a prescribed pressure 
boundary condition, while the right boundary was set as a gravity 


drainage condition to allow infiltrated water to exit the domain. 


The extension of the domain is long enough to allow the 
development of a perched aquifer in some of the simulated scenarios and 
small enough to allow for the use of fine grid discretization for the sought 
spatial scale and level of detail. The model domain is discretized into 


rectangular cells with variable horizontal spacing from 0.05 m near the 
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well, to 10 m at the outer boundary (r = 600 m). Vertically, the domain 
is discretized into layers of Az = 1 m, except for the layers that cross the 
well that had a finer discretization of Az = 0.5 m. Preliminary simulations 
showed that this level of refinement was sufficient to simulate changes in 
pressure and water table elevation near the well for the simulated 
scenarios, while maintaining the computational effort within reasonable 


levels. 


No flow boundary conditions are specified for the top and bottom 
boundaries of the domain and for the left boundary, which coincides with 
the center of the well or symmetry axis. The right boundary condition (r 
= 600 m) above the initial water table was specified as a seepage 
boundary (Lappala et a/., 1987), however the domain size is such that the 
perched aquifer that develops in some of the simulations never reaches 
this boundary. Below the initial water table, the boundary condition 
preserves the initial hydrostatic pressure distribution, so that water that 
reaches the boundary exits the domain. We consider a single injection 
well 57 m long and diameter D = 0.3 m with a screen (open section for 
flow) at the bottom of the well of length L = 13 m. Within the well section 
where water enters the domain, prescribed pressure increases with depth 
to simulate a quasi-hydrostatic pressure distribution along a 44 m water 
column stored inside the well casing. We simulated the operation of the 
system for two years: A one-year injection period, followed by a similar 
drainage period. Although the simulated scenarios are rather simple, they 
are a good approximation to real sites where the subsurface tends to 


exhibit quasi-horizontal stratifications. 
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However, the setup of the simulations does not consider potentially 
relevant processes such as: cyclic operation of the recharge system, air 
trapping and hysteresis due to application of wetting/drying cycles, 
heterogeneity within each hydrogeological unit (e.g. Sasidharan et al., 
2019; Sasidharan et al., 2020), capillary forces and potential temporal 
changes in properties due to clogging or chemical process (e.g. Glass et 
al., 2020). Adding some of those processes could potentially modify, at 
least in part, the results discussed below and, hence, should be 
considered for future work to improve on the simulation of the operation 
of real artificial recharge systems. Nevertheless, the results presented 


below represent a wide range of real drywell systems. 


Results 


Figure 4 shows 2D sections of simulated saturation, ¡.e. the fraction of the 
available pore space occupied by water, after 50 days and 1 year since 
injection started and after 1 year since injection stopped. The distribution 
of the simulated saturation for the homogenous case (S1) shows the 
development of a mounding zone directly below the injection well that 
extends almost 100 m from the well and remains almost constant during 
the injection period (first year). After injection stops, the water that was 
within this area drains under the action of gravity, so that by the end of 
the simulated period, saturation within this volume is smaller than 0.4. 
The results for scenario S2 are similar, but show a lower saturation value 
below the contact between HU1 and HU3 that is due to the lower K in HU3 
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and thus a lower vertical water flow. For the third scenario, that includes 
the presence of the low K unit (HU2) with K = 1.1 m/d, the lower 
infiltration rate in this unit produces as result the development of a 
perched aquifer that rapidly extends towards the outer boundary of the 
domain. The final extension of the perched aquifer is such that the total 
vertical infiltration is equal to the infiltration rate at the well. For example, 
considering a vertical velocity, V,, between the upper and lower aquifers 
of the order of the vertical hydraulic conductivity of the low K layer (0.01 
m/d) assuming that vertical flow away from the well is mainly 


gravitational, results in the following mass balance equation: 


where I is the infiltration rate and R is the maximum extension of the 
perched aquifer under steady-state conditions. For example, for S4, S5 
and S6 the infiltration rate is about 190 l/s (16550 m3/d); then R can be 


estimated as 725 meters, which is longer than the domain size. 
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Figure 4. Vertical 2D saturation profiles after 50, 365 and 730 days 
since injection started. Columns show results for the four simulated 


scenarios. 


The area below HU2 shows intermediate saturation values during 
and after the injection period, ¡.e. the soil in that area never gets fully 


saturated. At the end of the simulated period, saturation is still above 0.5 
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for all HU2. Figure 5 shows a 3D reconstruction of the results for this 
scenario, where the development of the perched aquifer can be easily 
identified after 5 and 365 days since injection started. Finally, the 
distribution of saturation for scenario S4, which considers a lower value 
of K for HU2 (K = 0.11 m/d), also shows the development of a perched 
aquifer during the injection period, which continues extending horizontally 
by the end of the simulated period due to the lower vertical infiltration 
rate that is allowed through HUZ2. 


S3 


5 days 365 days 730 days 


-100 , á -100 > E 
100 0 . 100 0 . 100 


100 200 Radius [m] Radius [m] 100 200 Radius [m] Radius [m] 100 200 Radius [m] 
Figure 5. Simulated saturation for scenario S3 after 5, 365 and 730 
days since injection started. 3D reconstruction from results of 2D 


axisymmetric simulations. 


Analyzing the hydraulic connection between the injection zone and 


the water table also provides a quantitative measure of the differences 
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between the scenarios. Figure 6 shows the simulated water pressure 
versus time for points located along the domain right below the water 
table (z = 140 m). The change in pressure is maximum for the 
homogeneous case where there ¡s a full hydraulic connection between the 
injection zone and the water table below. In this case, pressure increases 
by almost 10 meters below the injection well. The pressure increase is 
lower for the other scenarios because of the development of a perched 
aquifer or the presence of a partially saturated zone between the injection 
zone and the water table. The increase in pressure is the lowest for 
scenario S4, which considers the lowest K value for HU2. The increase in 
pressure for scenarios S2, S3 and S4 is due to the mounding of the water 
table produced by the increased vertical infiltration in the surroundings of 
the injection well, where the spatial variation of the piezometric head ¡is 
maximum. These differences can have practical implications for the 
design of monitoring systems that rely on the observation of variations of 
piezometric levels in monitoring boreholes and their interpretation (Gale 
et al., 2002). Relatively small changes in hydraulic conductivity values for 
some of the HUs can result in differences in pressure rise of up to a few 
meters during an initial period that might extend for a few months or 
years. In practice, small increases in pressure may be difficult to associate 
to the direct effects of artificial recharge; hence, they may not provide 
enough evidence of the effectiveness of the system and its benefits (Gale 
et al., 2002). 
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Figure 6. Pressure profiles along initial position of the water table (z = 


140 m). Colors indicate distance from injection well in meters. 


Figure 7 shows simulated saturation for scenarios S2, which 
considers only two HUs, and scenarios S5 and S6, which are similar to 
scenario S4, but with a low hydraulic conductivity layer that extends only 
up to 300 and 550 meters, respectively. The saturation front around the 


well is similar for all three scenarios at early time (5 days), but it is 
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significantly different at later times. Like in S4, the presence of the middle 
low K layer results in a perched aquifer that extends for a few hundred 
meters. In the case of S5, the perched aquifer extends up to the limit of 
the low K layer (300 meters), where water spills through a narrow vertical 
zone until reaching the water table. For scenario S6, the perched aquifer 
develops without reaching the limit of the low K layer (550 m) by the end 
of the injection period (365 days). Hence, the saturation distribution for 


S5 is similar to the one for S4, except for the section close to the border. 
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Figure 7. Simulated saturation for Scenarios S2, S5 and S6 after 5, 10, 
365 and 730 days since injection started. 
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Figure 8 shows the simulated piezometric head, h, distribution for 
Scenarios S2, S5 and S6. The contour for h=0 meters indicates the 
position of the boundary of the saturated zone that develops around the 
injection well and the original position of the water table. In all three 
scenarios, the upper saturated zone by the end of the injection period 
(365 days) extends only up to the position of the horizon between the 
upper aquifer and the underneath HU, lower K aquifer for S2, or low K 
middle unit for S5 and S6. 
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Figure 8. Simulate pressure head (m) for Scenarios S2, S5 and S6 after 


5, 10, 365 and 730 days since injection started. White line indicates 


where pressure head equal to O, ¡.e. boundary of saturated sediments. 


Figure 9 shows the normalized cumulative volume of infiltrated 


water that remains in the domain versus time (Table 3). For the first two 


scenarios, most of the injected water reaches the water table and leaves 


the domain quickly. The increase in storage produced by the injection is 


relatively minor in comparison to the total infiltrated volume after a short 


initial period. In a real case, this means that most of the infiltrated water 


can be recovered from the aquifer through the operation of deep wells 
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that reach below the water table. For scenarios S3 and S4, the proportion 
of the injected water that remains in the domain is higher, particularly for 
S4, because of the presence of the low K unit (HU2) that enhances the 
development of a perched aquifer and slows down the drainage of the 
infiltrated water after injection stops. This result can also have important 
practical implications since recovering water from a shallow perched 
aquifer can be more difficult than recovering it from the main existing 
aquifer (Bouwer, 1996; Bouwer, 2002). Alternatively, the implementation 
of a deep recovery system that extracts water directly form the aquifer 
should consider a potential long delay between recharge and the time the 


recharged water reaches the aquifer. 
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Figure 9. Normalized net recharged volume within model domain 


versus time. 
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Table 3. Mean infiltration rate and volume retained in domain at the 


end of the simulated period (2 years since injection started). Normalized 


values correspond to the ratio between the values and the maximum 


value in each column. 


Scenario 


s1 
S2 
S3 
S4 
S5 
S6 


Infiltration 
rate 


(m?/d) 


41 200 
17 320 
17 020 
16 450 
16 629 


16597 


Normalized 
Infiltration 
infiltration 
rate (l/s) 
rate 
477 1.00 
200 0.42 
197 0.41 
190 0.40 
192 0.40 
192 0.40 
Conclusions 


Retained 


Volume 2 years 


(x103 m3) 


106 
87 
558 
4 630 
2 086 


4 635 


Normalized 
retained 
volume 2 


years 
0.02 
0.02 
0,12 
1.00 
0,45 


1.00 


We performed detailed numerical simulations of an artificial recharge 


system composed of a single drywell that injects water into the vadose 


zone. 


The results of the simulation demonstrate that significant 


differences in the performance of the system in terms of injection rate, 


recharged volume and hydraulic connection between the injection well 


and the existing aquifer; can occur depending on the distribution of 
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materials that compose the subsurface and some of their characteristic 


hydraulic parameters. 


For practical purposes, it may be advantageous to implement 
injection systems that directly inject water into existing aquifers to avoid 
the potential development of perched aquifers, which can difficult the 
recovery of the injected water and delay the effective recharge to the 
main aquifer system. This may be particularly important in locations 
where the vadose zone includes a few low permeability units or the water 
table is deep as often happens in arid or semi-arid climates. Moreover, 
the design and implementation of monitoring systems to measure the 
performance of recharge systems can also be challenging in presence of 
heterogeneous hydrogeological units. In those cases, changes in 
piezometric head measured in monitoring boreholes might not provide 
relevant information to quantify recharge. Hence, potential cost savings 
due to using less expensive shallow infiltration systems, can be rendered 
of little worth due to the impossibility to recover the injected water within 
a reasonable time frame. These findings can be contradictory to 
recommended practices in some countries, where regulatory guidelines 
advice the injection of water into the vadose zone to profit from the 
potential filtration and natural chemical and biological mechanisms that 


can enhance the quality of the infiltrated water. 


Future work to extend the results of this study should consider other 
factors that can potentially impact the performance of a recharge system, 
such as clogging due to the transport of fine sediments, chemical 
reactions induced by the difference between the water quality of the 


injected and resident water, or the growth of microorganisms around the 
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injection point. Studies that analyze the potential performance of specific 
sites should consider a detailed characterization of the sediment 
properties, particularly unsaturated flow properties and anisotropy. 
Nevertheless, these additional improvements would affect only 
quantitatively the conclusions of this study but would not change their 


spirit. 
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Appendix: Unsaturated flow properties 


The unsaturated properties were set according to the van Genuchten 
model (van Genuchten € Nielsen, 1985), which includes the following 


equations: 


1 


e aa” (A1) 
k, = aer A (A2) 
D2 


265 
Tecnología y ciencias del agua, ISSN 2007-2422, 


Open Access bajo la licencia CC BY-NC-SA 4.0 15(1), 233-274. DOI: 10.24850/j-tyca-15-01-06 


(https: //creativecommons.org/licenses/by-nc-sa/4.0/) 


OPEN Que: 1) Check for updates 
Tecnología y 
CienciaszAgua 
With: 


C = la: h|é-=D 
D=1+ la: h)P 
Where: 


se: effective saturation or reduced water content (-), equal to s, = (0 — 
0,)/(n= 0,) 

k, : relative permeability (-). 

6: volumetric water content (-). 

9, : residual volumetric water content (-). 

n: porosity (-). 

h: pressure head (m). 

a: parameter of the model (1/m). 


f: parameter of the model (-). 


y =1-1/B 


In all simulations, we considered a loamy-sand soil with parameters 
summarized in Table Al. Figure A1l shows the corresponding water 


retention and relative permeability curves. 
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Table Ai. Parameters that define unsaturated flow properties for the 


loamy sand soil considered in all simulations (Lappala et a/., 1987). 


Parameter Value 
Saturated Hydraulic Conductivity EEG 
K (m/d) 
Porosity, n (-) 0.41 
a (1/m) 12,4 
B 2.28 
Y 0.56 
Residual Water Content, O, 0.057 
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Figure Al. Water retention and relative permeability curves for the 
loamy sand soil considered in all simulations. 
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